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Investigating the Andromeda Stream: I. Simple Analytic 
Bulge-Disk-Halo Model for M31 



O 
O 

(N 

Sh 

Oh 
<C 
O 



(N 
> 
O 

(N 

o 
in 
o 

6 



13 



J. J. Geehan 1 *, M. A. Fardal 1 ' 2 , A. Babul 1 , P. Guhathakurta 3 

1 Dept. of Physics & Astronomy, University of Victoria, Elliott Building, 3800 Finnerty Rd., Victoria, BC V8P 1A1, Canada 
2 Dept. of Astronomy, University of Massachusetts, Amherst, MA 01003, USA 

3 UCO/Lick Observatory, Dept. of Astronomy & Astrophysics, Univ. of California, 1156 High St., Santa Cruz, CA 95064, USA 



Accepted 2005 November 10. Received 2005 October 30; in original form 2005 January 13 



ABSTRACT 

This paper is the first in a series which studies interactions between M31 and its 
satellites, including the origin of the giant southern stream. We construct accurate yet 
simple analytic models for the potential of the M31 galaxy to provide an easy basis for 
calculation of orbits in M31's halo. We use an NFW dark halo, an exponential disk, 
a Hcrnquist bulge, and a central black hole point mass to describe the galaxy poten- 
tial. We constrain the parameters of these functions by comparing to existing surface 
brightness, velocity dispersion, and rotation curve measurements of M31. Our descrip- 
tion provides a good fit to the observations, and agrees well with more sophisticated 
modeling of M31. While in many respects the parameter set is well constrained, there is 
substantial uncertainty in the outer halo potential and a near-degeneracy between the 
disk and halo components, producing a large, nearly two-dimensional allowed region 
in parameter space. We limit the allowed region using theoretical expectations for the 
halo concentration, baryonic content, and stellar M/L ratio, finding a smaller region 
where the parameters are physically plausible. Our proposed mass model for M31 has 
Mbuige = 3.2 x 10 10 M , M disk = 7.2 x 10 10 M©, and M 200 = 7.1 x 10 11 M , with un- 
corrected (for internal and foreground extinction) mass-to-light ratios of M/Lr = 3.9 
and 3.3 for the bulge and disk, respectively. We present some illustrative test particle 
orbits for the progenitor of the stellar stream in our galaxy potential, highlighting the 
effects of the remaining uncertainty in the disk and halo masses. 
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1 INTRODUCTION 

The detritus resulting from the disruption and assimilation 
of a satellite galaxy by its larger host is a veritable treasure 
trove of clues and information about the galaxy assembly 
process, the nature of the merging subunits, as well as 
the dynamical properties of the parent system. Typically 
in the form of stellar streams, this detritus can remain 
spatially and kinematically coherent in the halo of the 
larger galaxy fo r sever al billion years iJohnston et aljll996t 
iHelmi fc White! I1999T) . Identifying and quantifying these 
features in detail, however, is a challenging task as it often 
entails resolving and measuring properties of individual 
stars in the stellar halo of the galaxies. At present, such 
studies are only feasible for the galaxies in the Local Group. 
The Milky Way stellar halo has been the focus of detailed 
scrutiny for many years now resulting in the detection 
of a number of coherent features in the star counts (see 
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Yannv et alJ l2003t iMaiewski et al 



l20(5lLaw^olin^ton^^M^ewsk1f 

therein). Unfortunately, as noted by 
our vantage point within the Galaxy 
cates the interpretation of these features. C onsequently, 



200 
2005 



iMartin et alJ 

and references 



Vewberg et all J2002 
greatly compli- 



the recent uncovering 


>; of rich substructure dlbata et alJ 


20011 IFereuson et al. 


l2002t iMcConnachie et alJ 120031: 


Morrison et all .2003: 


iMerrett et alJ 120031: IFereuson et alJ 


2004 iGuhathakurta et alJ 12004 IZucker et alJ 12004 


Merrett et al J l2004h in the stellar halo of the Andromeda 



galaxy (M31) is highly tantalizing. 

The most striking of these features is the giant southern 
stellar stream extending out from the south-ea stern part of 
M31's disk. First reported bv llbata et al] i200ll) . the stream 
stars have since been targeted for careful photometric and 
spectroscopic analyses resulting in the determination of the 
distance to the stream at various locations along its length 
iMcConnachie et alj[2003l) as well as the measurement of 
stellar kinematics in a num ber of locations dlbata et all200j : 
IGuhathakurta et al]l200'3) . Jointly, these observations indi- 
cate that the giant southern stream extends away from us 
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below the disk of M31, out to a distance of ~ 100 kpc away 
from the M31 center. There is also a strong velocity gradient 
along the body of the stream, with the outer regions being 
at rest with respect to M31 while the inner regions approach 
us at 300 km s _1 with respect to M31. 

Intense observational effort notwithstanding, the source 
of the giant southern stream has yet to be identified. We 
do not know whether the progenitor survives and if so, 
where it is, although several faint features in the inner 
halo of M31, such as t he newly detected satellite And VIII 
jMorrison et alJ 120031) . have been touted as possibilities. 
One approach to identifying the progenitor of the stream, 
or at least its fate, is to reconstruct its orbit using the 
stre am propert i es as constraints. Both llbata et alJ J2004T) 
and iFont et "all i2004T) have initiated efforts in this direc- 
tion. As we shall illustrate, both in this and in Paper II 
jFardal et al.ll2005l) . a detailed study of the stream dynam- 
ics and of its progenitor orbit, however, requires a realistic, 
reasonably accurate model for the potential of M31 that is 
also preferably easy to use and straightforward to alter for 
experimentation purposes. 

The simplest approach to describing the potential 
of M31 is by modeling the mass distribution in the 
galaxy. The earliest efforts at doing so date back to 
iBabcockl il938. 1938ft. fo l lowed bv_ a ser i es of i nvestiga- 
tions b y IWvse & Mavall il942h . iKuzminl dl943l). iKuzmirJ 
lll952h . ISchwarzschildl Jl954h . and ISchmidtl lll957l) over 
the course of the next two decades. Subsequent im- 
provement in the quality of photometric and spectro- 
scopic data in the 70s and the 80s led to a r e assess - 
ment of the mass mod el s bv iDeh arveng fc Pellet) l W7t ) , 
iMonnet fc Simienl lll977MSimien. Pellet, fc Monnetl Jl979f) . 
and iKentl 1119891) . As heroic as these efforts were, these 
models were constructed at a time when the shape and 
the parameters describing dark halos were largely un- 
known and unconstra ined. Most recently, the problem 
has been revisited bv |Klvpin 1 Zhao. fe Somerville | l|2002|). 
IWidrow. Perrett. fc Suvul J2003D and IWidrow fc Dubinskil 
i2005l) . However, the mass models put forth by these 
groups, though sophisticated, are neither simple nor easy 
to use especially for the purpose s of orbit calculations. The 
IWidrow. Perrett. fc Suvul (120031) mass models, for example, 
are specified in terms of a set of distribution functions. While 
these distribution functions can be used to compute the mass 
density distribution for M31 as well as the associated grav- 
itational potential, the derivation is implicit and must be 
solved for iteratively. In other words, a closed-form analytic 
description of the density and the potential is not available. 
Similarly, the Klypin et al. model is also not available in 
closed analytic form. 

In this paper, we attempt to remedy this by presenting 
a reasonably accurate yet simple analytic description of the 
M31 mass distribution. Our principal aim is to arrive at a de- 
scription of the gravitational potential that is both suitable 
for the purposes of computing satellite orbits and sufficiently 
transparent that it can be easily altered for the sake of ex- 
perimentation. We achieve this by decomposing the mass 
distribution in M31 into four components: the central black 
hole (BH), the bulge, the disk and the extended halo. The 
components are modelled using well-known functional forms 
whose associated potentials are easy to compute. The appro- 
priate values of the associated structural parameters are de- 
termined by requiring that the model is in agreement with 



the observed M31 rotation, surface brightness, and bulge 
velocity dispersion profiles. In § |5J we describe the analytic 
functional forms we have chosen for the four components. In 
§ [3] we review the M31 observations used to constrain the 
structural parameters, and visually compare to our best-fit 
model. In § 3] we examine the allowed region of parameter 
space, which shows both strong constraints in some direc- 
tions and near-degeneracies in others. We discuss additional 
physical constraints on the solutions besides those used in 
the fit. Two "constrained best-fit" solutions resulting from 
this discussion are presented along with those of our formal 
best-fit solution in Tabled we direct readers interested only 
in the final product to this table. In §|KJ we briefly illustrate 
the relevance of our present efforts, by computing sample 
orbits for the progenitor of the giant southern stream. We 
discuss the effect of the uncertainty in the potentials on the 
orbits, and compare our orbits to those obtained in some 
simple analytic potentials that have been recently used for 
this purpose. A summary of our key results is presented in 
§ H2 For completeness, we note that our adopted cosmol- 
ogy is a spatially flat ACDM universe with Sl m = 0.14/i -2 , 
Sit = 0.024/i~ 2 , and a Hubble constant of h = 0.71 in units 
of 100 km/s/Mpc. 



2 COMPONENTS OF OUR M31 MASS MODEL 

In principle, a detailed mass model of M31 ought to comprise 
of several components: the central black hole, the nucleus, 
the bulge, the bar, the spheroid, the thick and thin stellar 
disks and the thin gaseous disk. In the interest of simplicity, 
we restrict ourselves to only four components: the central 
BH, the bulge, the disk and the extended halo. Moreover, 
our use of simple analytic functions to describe these com- 
ponents necessitates making some simplifying assumptions. 
We discuss these as well as our choice of functional forms 
below. 



2.1 The Central Black Hole 

The very central region of M31 is comprised of a central BH 
and a distinct small-scale stellar nuclear component that 
is photometrically and dynamically separate from the bulge 
and t he large-scale galactic disk (see lKormendv fc Richstond 
Il995l and references therein). The most re cent estimate of 
the B H mass is M. = (5.6±0.7) x 10 7 M dSalow fc Statlerl 
2004), while the mass of the stellar com ponent is estimated 
to be < 2xl0 7 M fl , within ~ 10 pc feacon et al.ll200ll : 
|Pensll2002l:ISalow fc Statleill2004l) . 

In terms of the total M31 potential, the two nuclear 
components are only relevant on scales r < 20 pc and for 
all practical purposes, can be ignored in the context of the 
calculations of interest to us. However, to facilate wider use 
of our M31 model, we include the dynamically more impor- 
tant of the two, the BH component, in the form of a point 
mass located at the center of our galaxy, and fix its mass 
to the value quoted above. The inclusion of this component 
is necessary, for example, in order to correctly model the 
velocity components in the nucleus of M31 (cf. Figure [U. 
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2.2 The Galactic Bulge 

We model the bulge of M31 as a spherically sym metric mass 
distri bution represented by a Hernquist profile jHernauisd 
1990): 



Pb(r) = 



M b 
2-Kr'l 



{r/r b ){l + r/r b f 



(1) 



where Mb is the total mass of the bulge and rb is its scale 
radius. The mass profile and the potential corresponding to 
this density distribution are 



M b (r) 



Mr) = - 



M b r 2 
(rb + r) 2 ' 
GM b 
r b + r' 



(2) 
(3) 



We con s idered using the more general density profile of 
iDehnenl l|l993f) where p b (r) oc (r/r 6 )~ 7 (l + r/r 6 ) 7-3 but 
found that the results for 0.1 < 7 < 2 were all equally good. 
For simplicity, we have opted to use the 7=1 Hernquist 
profile. We also note that the both the bulge mass and the 
bulge mass-to-light ratio are quite insensitive to variations 
in 7 within the range noted. This gives us confidence that 
the bulge parameters that we quote in Table |5J and more 
generally our bulge-disk decomposition, are robust regard- 
less of the precise nature of the inner density structure in 
the bulge. 

On the other hand, using a spherical model for 
the bulge is clearly a simplification. Detailed modeling 
based on surface photometry suggests that at the very 
least, the bulge ought to be mode lled as an oblate 
spheroid w i th axis ratio of ~ 0.8 ||Kent| Il983t Il989t 
iPend 120021 IWidrow. Perrett. fc Suvl 12003ft and that in 
actuality, it is almost certainly triaxial jLiridloladJ E256t 
Starkll977tlKentll989l:IStark fc BinneJl994tlBermar]l20oJ: 



Berman fc Loinardl 12002'). In opting to avoid the compli- 
cations of using a non-spherical mass distribution, we are 
primarily guided by the purpose of this paper, which is to 
construct orbits to model the giant southern stream. This 
stream has a pericenter estim ated to be at ~ 2 to 4 kpc, 
and an apocenter of ~ 100 kpc Jlbata et all2004l : lFont" et alJ 



2004), so it spends little to no time in a region where the 
asphericity of the bulge could have any dynamical effect. 
After all, the equipotential surfaces tend to spheres at large 
radii even if the density distribution does not. However, if 
our models are used to treat dynamics at smaller radii, the 
reader should keep this simplification in mind. 



2.3 The Galactic Disk 



In order to model 
gin by assuming, as 



IWalterbos fc KennicutB (U98 
bution can be described by an exponential surface density 
profile: 



the M31 galactic disk, we be- 
suggested by the observations of 
that the disk mass distri- 



E d (7?) = T, oe - R/Rd 



(4) 



where Eo is the central surface density, R d is the scale length 
of the disk, and R is the distance from the centre of M31 
in the plane of the disk. This corresponds to a disk mass 
within a sphere of radius r of 



M d {r) = 27rE 0j Rl [l - (1 + r/R d )e- r/R ^ 



(•>) 



Here, r = \/ R 2 + z 2 is the distance from the centre of M31 
and z is the distance perpendicular to t he plane of the disk. 

F or an infinitesimally thin disk, iBinnev fc Tremaind 
lll987l) give the expression for the potential as: 



$d(R,z) 



-2nGY, Ri 



J (kR)e- k ^ dk 
[1 + (kR d ) 2 } 



213/2 ' 



(6) 



and in the plane of the disk, this thin disk potential implies a 
circular velocity profile given by l|Binnev fc Tremainei 198^ 1 : 



V c 2 d (R) = 4nGZ R d y 2 [Io(y)K (y) - h(y)Ki(y)] , (7) 

where y = R/2R d and (Jo, Io, Kq, ii, Ki) are Bessel func- 
tions. 

Below, we shall also have occasion to use an extremely 
simple toy model that we will refer to as the "spherical disk" 
model. In this model, we assume that the disk mass within a 
sphere of radius r is distributed in a spherically symmetric 
fashion rather than concentrated in a thin disk. Then the 
corresponding gravitational potential is simply 



1 - e" 



r/Rd 



(8) 



Unless otherwise stated, the results presented in this paper 
will be based on the thin disk model. 



2.4 The Extended Dark Halo 

Finally, we consider the extended dark matter halo of M31. 
We assume that this component can be adequately modeled 
as a spherically symmetric system. We recognize that the 
issue concerning the sphericity of the halo is a contentious 
one, but given that it remains observationally unresolved 
even in the case of the Milky Way, we feel it more in keeping 
with the spirit of our approach to adopt the simplest model. 
To describe t he run of density with radius, we adopt the 
NFW profile (iNavarro. Frenk. fc Whitelll99ffl : 



Ph(r) 



S c p c 



{r/r h ){l+r/r h y 



(9) 



where p c = 277.72ft 2 Mg/ kpc 2 is the present-day critical 
density, h = 0.71 is the Hubble constant in units of 100 
km/s/Mpc, S c is a dimensionless density parameter, and 
is the halo scale radius. N-body simulations based on the hi- 
erarchical clustering scenario for structure formation within 
the cold dark matter cosmogony suggest that spherically av- 
eraged halo density profiles are well described by the above 
profile. There is much debate over the exact exponent of 
the density in the inner cusp, but in our specific instance, 
this will not matter since the potential within the region 
in question will be dominated by the disk and/or the bulge 
components. 

The mass profile and the corresponding potential for an 
NFW halo are given by: 



M h (r) =4-KS c p c ri 



hi 



r + r h 



r + r h 






r + r h 




Th 



(10) 

(11) 
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3 SPECIFYING THE M31 MODEL 
PARAMETERS 

The BH-Bulge-Disk-Halo mass model described above has 
a total of 6 structural parameters remaining to be specified: 
8 C , Thi Mb, rt, So, and Rd, or in other words a scale radius 
and normalization for each of the bulge, disk, and halo com- 
ponents. We constrain these parameters using a number of 
extant observations of M31. 

Before discussing the data and our efforts to construct 
the mass model, we briefly review what is known about the 
configuration of M31. The Andromeda galaxy lie s at a dis- 
tance of 784±24 kpc JStanek fc Garnavicblll998ft from the 
Milky Way, and has a mean rad ial velocity of —300 ± 4 
km/s ide Vaucouleurs et al. lll99lT) . We look at the galaxy 
from below, as is evident from images of the disk dust 
lanes projected onto the bulge. We take the inclination to 
be 77°. This is the generally ac cepted value, although es- 
timates range from 74° to 7 9 ° (iRubin fc D'Odoricojll969t 
IWalterbos fc Kennicutdll988t iMa. Peng, fc Gulll997TI . The 
galaxy spins counterclockwise on the sky. We take the po- 
sition angle of the disk to be 37°. The position angle of 
the bulge major axis app ears to be offset from that of the 
disk by ~ +10° to +20° jMcElrovHl983 lKormendvlll988D 
and additionally, the bulge isophotes have distinctly "box- 
shaped" appearances with ellipticities that increase with 
radius. These latter features, as well as the anomalously 
high non-circular gas velocities in the inner disk, have gen- 



erally been interpreted as si 


gnatures of 


jLindbladl Il956 


: Starkl 119771 


iKend Il98£ 


1994: iBermanl 


2001: iBerman & Loinarc 



In addi- 
tion, analyses of both the HI and the light distribution 
in M31 disk suggest that the disk is significantly warped, 



especially in the outer reg i ons fcf. ISawa fc Sofud 
|lnmmen^£amrjej^j^iy^}82^^ 
Il988l: [ 



Il988l: iMorris. Reid. Griffiths fc Pennvlll99l . iBraunl dl99lTl 

argues that in addition to this warping, the inclination of 
the M31 HI disk also varies as a function of both radius and 
azimuth and that it may even be flaring toward the outer 
regions. The latter is not altogether surprising in light of 
the number of satellites that appear to be interacting with 
the galaxy; the flaring of the disk is one of the generic fea- 
tures of such interactions (cf. Fardal, Babul & Barnes, in 
preparation). 

Over the past two decades, detailed studies of M31 have 
generated a large volume of kinematic and photometric data. 
For practical considerations, we restrict ourselves to a lim- 
ited number of data sets. Specifically, we constrain our mass 
model using both the M31 major-axis and average surface 
brightness profiles, the bulge major and minor axis velocity 
dispersion profiles, the disk rotation profile, and the dynam- 
ically derived total mass estimates in the intermediate and 
outer halo. These are described more fully in § I3.1H3.4I In 
order to relate our mass model to the light profiles, we are 
required to introduce two additional parameters: (M/Ln) b 
and (M/Lr)j, or the mass-to-light ratios of the bulge and 
the disk. We assume that these two ratios are constant over 
the entire bulge and disk; this is probably not exactly cor- 
rect but should b e reasonably valid, given the modest color 
gradient in M31 (IWalterbos fc Kennicutdll987l) . 

The inclusion of these mass-to-light ratios increases the 
number of parameters to 8. We tackle the problem of fit- 
ting these parameters in two steps. First, we fit the sur- 



face brightness profile to obtain the scale radii r& and Rd 
of the bulge and disk, and the luminosity normalization of 
these components. The abundance and quality of the sur- 
face brightness data is such that it would dominate the fit 
of the scale radii in any case. Also, the x 2 from this fit- 
ting step is poor, which merely tells us that we have not 
managed to fit all the details of M31's visual structure; in 
absolute terms, the size of the discrepancies is actually fairly 
small. Then, we incorporate the dynamical information to 
fit the four remaining parameters (the bulge and disk mass 
normalizations Mb and So, and the two halo parameters Vh 
and (5 C ). In contrast, the % 2 from this fitting step is actually 
meaningful, and we use it to constrain the allowed region of 
parameter space. 

We present the parameters for our best-fit solution in 
Table |21 and as we introduce below the set of observations 
that we use to constrain the model, we also describe how 
well our best-fit model fares in comparison. We defer until 
§ 2] any discussion of errors, correlations, and possible de- 
generacies between the parameters, and the correspondence 
between the allowed terrain of parameter space and physical 
properties of M31. 



3.1 M31 Surface Brightness Data 

Given its proximity, M31 is one of the few galaxies visible 
to the naked eyes. Its existence was first documented by the 
Persian astronomer Abd al-Rahman Al-Sufi in his treatise 
on stellar astronomy titled "Kitab al-Kawatib al-Thabit al- 
Musawwar" (Book on the Constellations of the Fixed Stars), 
published in AD 964, where he both identified its position 
in the sky and summarized his observations. Though sub- 
1982; j ec ^- Q £ some speculation and study over muc h of the ensuing 
1987, millemium, it was Hubble's epochal study (iHubble Ill929l) 
establishing the true nature of M31 as an extragalactic stel- 
lar system (galaxy) that underlies the considerable effort 
devoted to the study of M31 since the 1930s. As a result, 
there is no dearth of high resolution photometric data on 
the system. Here, we will use both the azimuthally-averaged 
(or "global") surface brightness profile and the major-axis 
profile to constrain our models. The global profile has the 
advantage of minimizing effects from localized structural fea- 
tures, like the spiral arms, which can introduce bumps and 
wiggles in the light curve along any one direction. In the 
left p anel of Fig ure ^ w e plot t he J?- band global light pro- 
file of IWalterbos fc Kennicutil Jl987l) . We use the .R-band 
profile because we expect it to be less affected by dust ex- 
tinction and stellar population variations than bluer bands. 
For reasons discussed below, w e find it necessary to shift the 
IWalterbos fc Kennicutd dl987l) major axis surface brightness 
profile by —0.1 mag, and since the major axis and the present 
global light profiles are derived from the same photographic 
plates, we apply the same shift to the global profile as well. 
Walterbos & Kennicutt compute the mean surface bright- 
ness at radial coordinate r by averaging the light in an el- 
liptical annulus formed by taking a face-on circular ring of 
radius r and thickness 8r, rotating it to an inclination of 
77°, projecting the resulting structure onto the plane of the 
sky, and aligning the major axis of the resulting elliptical 
annulus with the major axis of M31. When we compute the 
global surface brightness for our model, we follow the same 
procedure. 
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Figu re 1. Left panel: The global surface brightness profile of M31 in log-linear projection. The red squares represent the data of 
Walt erbos fc Kennicut j ll987f) . corrected as described in the text. The blue lines show the results of our best-fit model. The dotted line 
is the bulge contribution, the dashed line the disk^and the thick solid bl ue line is the total. For comparison purposes, we juxtapose the 
observed major axis surface brightness data of Walterbos & Kcnnicutt ( 1987) from the right panel as the thin red line. Right panel: 
The symbols trace the composite, observed M31 major axis surface brightness profile in log-log projection. For comparison, we plot the 
observed global light from the left panel as a thin red line. As in the left panel, the blue lines represent the major axis light profile for 
our best-fit model. 



The IWalterbos fc Kennicutl lll987l) global light pro- 
file only extends down to just under 1 kpc along the ma- 
jor axis, and therefore just begins to probe the region 
where the light is dominated by the bulge. On the other 
hand, there have been a number of studies of the light 
profile along the major axis that probe all the way in 
to the very central nuclear region. We, therefore, also in- 
clude the major-axis profile among our set of constraints. 
We plot this profile in the right panel of Figure ^ This 
combi nes the nuclear region HST data f rom l Lauer et alJ 
(1993), interme diate radii data fromlKentl jl983i) . and large 
radii data from IWalterbos fc KennrcuTTll98^K To convert 
Kent's r-band data and Lauer et al.'s l/-band data into 
an 7?-band light p rofile, we use colours V — R ~ 0.7 5 
jTenies et aU 11994) and r - R « 0.35 iJorgensenlfl994l) . 
As a check, we have tested these color transformations by 
running a stellar bu rst model in the PEGASE.2 popula- 
tion synthesis code dFioc fc Rocca-Volmerarigel Tl997). For 
a burst time of 12 Gyr, we find agreement within ±0.1 mag 
between our adopted offsets between V, r and R, as well 
as other measured colours of M31's bulge presented in Ta- 
ble Q Additio nally, we found it necessary to translate the 
IWalterbos fc Kennicutd lll987t) data by —0.1 mag in order 
to bring it in line with the other two data sets. This level of 
fine-tuning is not surprising given an uncertainty of ~ 0.1 
mag in the colours as well as ~ 0.1 mag zero-point errors in 
the surface brightness measurements. 

The major-axis light profile strongly indicates that the 
stellar nucleus (r < 0.01 k pc) and the bulge are two sep- 
arate components (see also iKormendv fc Bende"rlll999F) . As 



noted previously, we do not attempt to model the nucleus 
because its mass is too small to affect the masses of our 
galaxy components or orbits in the halo region. The major- 
axis profile also has a series of prominent bumps not present 
in the global light curve, due to dust lanes, disk warps, and 
spiral arms that intersect the major axis. Azimuthal aver- 
aging minimizes the presence of such features in the global 
profile. 

We compute the model major-axis surface brightness 
profile as fi(r) = /it(0.9r) + jj,d(r), where the disk and bulge 
surface brightness profiles in the ii-band are given by 

ix d {r) = S oe - r/ ^se C (i), (12) 

M»-) = 2(£) 1 jH p b {x)-^=dx. (13) 

Notice that in computing the surface brightness at projected 
distance r along the major axis, we sum the disk surface 
brightness at that projected distance with the bulge sur- 
face brightness at 0.9r. We do so in order to compensate 
for the fact that in our model the light distribution is cir- 
cularly symmetric while the observed light distribution is 
ellipsoidal. The factor 0.9 comes from our requiring that the 
area enclosed by an elliptical isophote at distance r along the 
major axis is the same as the area of our equivalent circle. 
This "correction" results in our circularly symmetric bulge 
light profile being stretched outward. 

When fitting the major axis and the azimuthally- 
averaged surface brightness data, we assign an uncertainty 
of ±0.1 mag arcsec -2 to each data point. The random er- 
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Table 1. Observed M31 Bulge Colors 



Bulge Color Reference 



r - Rzz 0.35 CJorgensen 1994) 

B - r S3 1.28 (Kent 19871 

V - R ss 0.75 fTenies et al. 19941 
B-Ra 1.7 fWalterbos fc Kennicutt 19871 

B - V S3 0.97 CWalterbos fe Kennicutt 1987: Tenies et al. 19941 

V - K S3 3.38 CPritchet 19771 

V - I fa 1.34 fLauer et al. 19931 

g - r S3 0.55 (Hoessel fc Melnick 19801 

v - g S3 0.78 (Hoessel fc Melnick 19801 

* All colors have an uncertainty of ±0.1 

rors are probably a strong function of radius since they drop 
with surface brightness, but there are significant systematic 
errors from zero-point shift s and color correction terms in 
all parts of these diagrams dWalterbos fc Kennicufrjll987ll . 

Furthermore, there are many features that cannot be 
fit using our simple model, and some of these are primarily 
features in the luminosity (e.g., dust lanes) and not in the 
mass. Our goal is simply to get a close approximation to the 
surface brightness at all radii, not to follow every feature. To 
this end, we use a standard \ 2 fitting routine to obtain the 
best-fit parameters, but ignore the actual value of x 2 which 
indicates a poor fit (reduced \ 2 ls about 7). 

The two panels of Figure show the global and major 
axis light profiles for our best-fit mass model for M31. Over 
most of the region of interest, the model agrees with the data 
within 0.2 mag. At first glance, the systematic offset of ~ 0.2 
mag at r > 10 kpc between the model and the data global 
surface brightness profiles may seem disconcerting. This off- 
set is due to the fact that for r > 10 kpc, the observed 
global light profile is brighter than the major-axis profile by 
as m uch as 0.5 mag at r — 20 kpc dWalterbos fc Kennicutt) 
Il987f) . To facilitate comparison, we show the trace of the 
observed major-axis and the global light profiles in the left 
and right panels of Figure respectively. If the disk was in- 
deed an infinitesimally thin, perfectly axisymmetric system 
and the bulge truly spherically symmetric, the major-axis 
and the global light profiles of the disk would be indistin- 
guishable. The observed differences in the light profiles are 
due to the increased contribution to the global light profile 
by the bulge due to its ellipsoidal shape, the warpage in the 
disk, and its finite thickness. Since our mass model does not 
include these additional features, we do not expect to be 
able to model the offset. Instead, our best-fit model "splits 
the difference" and settles in between the two. 



and thick dot-dashed curves, respectively). Our model re- 
sults are in excellent agreement with the observations and 
we remind the reader that the corresponding model param- 
eters are listed in Table [5] 

In the top left panel, we also show the results of a model 
that we have dubbed the "spherical disk" model. The halo 
and the bulge components in this model are exactly the same 
as in the best-fit axisymmetric model. In fact, even the total 
radial mass profile for the disk component, Md(< r) is the 
same in the two cases. In the spherical disk case, however, 
this mass is assumed to be distributed in a spherically sym- 
metric fashion, as opposed to planar and axisymmetric. As 
Figure |21 illustrates, the rotation profiles for the two models 
(and hence, the potential, at least in the plane of the disk 
of the axisymmetric case) are similar with minor differences 
being solely due to the different geometries of the disk mass 
distribution. Since orbit calculations in purely spherical po- 
tentials are much faster to compute, we use the spherical 
disk model to narrow the space of potentially acceptable 
satellite orbits and then, use the axisymmetric model to ex- 
plore this limited space further. This latter step is necessary 
because as we demonstrate in H4.2I the nature of the disk 
potential, whether axisymmetric or spherically symmetric, 
impacts the satellite orbits in non-trivial ways. 

In the bottom left panel, we juxtapose our best-fit ax- 
isymmetric model against the results for a model where 
the disk mass-to-light ratio has been constrained to be 
M/Lr = 3.3, which is lower than the mass-to-light ratio 
of the disk in the best-fit case. In the top right panel, we 
juxtapose our best fit-axisymmetric model against the re- 
sults for a model, which we will refer to as the "maximum 
disk" model, where the disk mass-to-light ratio has been 
constrained to be M/Lr = 6.3, a value that is much higher 
than the M/L of the disk in the best-fit case. The parame- 
ters for the two constrained models are also given in Tabled 
and reasons for considering these two additional models are 
discussed in iJ4.ll Here, we merely note that total rotation 
profiles for the formal best-fit model, the constrained-M/L 
model and the maximum disk model are virtually indistin- 
guishable. Fixing the disk mass-to-light ratio forces the disk 
mass higher or lower relative to the best-fit model. The cor- 
responding change in the disk contribution to the total ro- 
tation profile is entirely compensated by an opposite change 
in the halo contribution. This playoff between the disk and 
the halo is an example of the degeneracy that, in the ab- 
sence of any additional constraints, is a fundamental source 
of uncertainty affecting all efforts to construct mass models 
of disk galaxies. We discuss the nature of, as well as our 
efforts to break this degeneracy in H4.ll 



3.2 Disk Rotation Curve 

Figure [5] shows the rotation profile of the M31 disk. As 
with the photometric data, there is no dearth of optical 
and radio observations of the rotation data for the M31 
disk. We use the smoothed, co mposite rotation curve of 
IWidrow. Perrett. fc Suvul i2003ll . which was b ased on the 
HII emission regi on observatio ns of lKend jl989f) and the HI 
measurements of iBraur] dl99lll . 

The two left panels and the top right panel of Figure [5] 
also show the total disk rotation curve of our best-fit model 
(thick solid curve), along with the contributions from the 
bulge, disk and halo components (thick dotted, thick dashed 



3.3 Bulge Velocity Dispersion 

In Figure [3] we show the projected bulge velocity disper- 
sion profile along the bulge major (PA=45° — 55°) and minor 
(PA=135° — 160°) axes. The former is plotted in the upper 
panel while the latter is shown in the lower panel. The filled 
data points (circles) with error bars show the observed dis- 
persion that we are using in our fitting procedure. These 
correspond to the measurements bv lKormen d3 lll98dl along 
PA=55° (major) and PA=160° (minor). We only consider 
measurements at r ^ 0.02 kpc to avoid complications due to 
the BH. 
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Figure 2. The top left panel shows our best-fit axisymmetric (thick red lines) and spherical disk (thin cyan lines) model rotation 
curves, in comparison to the M31 disk rotation profile. The solid curves show the total model rotation curves. The short dashed lines 
show the contribution to the total from the disk, the dotted lines, from the bulge, and the dot-short dashed lines, from the halo. The 
spike in the curves at r = is due to the presence of the BH. The bottom left panel compares the best-fit axisymmetric and constrained 
(M/L=3.3) disk models. As in the top panel, the thick red lines correspond to the best-fit axisymmetric model. The thin blue lines 
show the results for the M/L=3.3 model. The line types are the same as in the top panel. The top right panel compares the best-fit 
axisymmetric (thick red lines) and the maximum disk (thin magenta lines) models. Line types are the same as discussed above. The 
bottom right panel compares our best- fit model rotation curve (thick red lines) against those of recently used models: the lBekki et alj 
l200l|) model used by Font et al. (cyan long-dashed line) , as well as the high (green short dashed line) and the low (green dot-dashed 
line') iKlypin. Zhao. &: Somerviliel 120021) models used by Ibata et al. 



For comparative purposes only, w e also include on the 
plot data points from lKormendvl Bl988fl taken along PA= 38° 



(open circles) as well as the data from lKormendv fc Bender! 
( 1999) (open triangles). The latter extends the velocity dis- 
persion profiles well into the nucleus of M31 and the trend 
in the measurements towards smaller radii reflects the in- 
creasing gravity due to the BH. Focusing on the PA=38° 
data from lKormendvl l|l988i) (open circles) , we note that this 
latter set is indistinguishable from that for PA=55° in the 
region of overlap, suggesting that small variations in the po- 
sition angles along which the measurements are made or the 
fact that these measurements may not be mapping out the 
velocity dispersion profile precisely along the "true" bulge 
major axis does not appear to introduce any significant bi- 
ases. 

For completeness, we also show the data points 
(open squares) corresponding t o the smoothed profile of 
IWidrow. Perrett. fc Suvul J2003t) . which itself is based on 
data fro'rn TMcl^ho^ " (^?83) . We only show the data only out 
to 1 kpc, since beyond this point we are concerned about 
the dynamical effect of bulge rotation and/or contamina- 
tion from the disk. There is one aspect of this data set that 
is slightly disc oncerting. In the r egion of overlap along the 
major axis, the lKormendvl dl988l) data i mplie s ~15% higher 
velocity dispersion than the iMcElrov! (119831) observations. 
We do not understand the origin of this disagreement. 

The solid curves in the two panels show the velocity 



dispersion for our best-fit mass model. The results match 
fairly well to the data points actually used to constrain our 
model. In keeping with the observations, the model profiles 
decline towards increasing radii beyond r = 0.2 kpc thoug h 
the decline is not as steep as indicated bv iMcElrovl (1983). 
Moreover, although the mod el was not constrained to fit the 
iKormendv fc Benderl fll999f) data, it nonetheless reproduces 
the gentle decline towards the center inward of 0.1 kpc, the 
trough at ~0.01 kpc, and the subsequent rise into the center. 

In computing the model velocity dispersions, we have 
made three additional simplifications, in keeping with our 
basic assumption of a spherically symmetric bulge. These 
additional assumptions are: (a) The bulge velocity disper- 
sion is isotropic; (b) the disk potential in the central 1.2 kpc 
can be approximated by that of the spherical disk; and (c) 
the bulge rotation can be neglected. Then the true (unpro- 
jected) bulge velocity dispersion is 



Pb(r' 



d$>tot , / 

~^rr dr ■ 

or 



(14) 



Pb(r) 

The observations, however, yield the luminosity-weighted 
pr ojected velocity dispersion. The equivalent model profile 



is (jsimien. Pellet, fc Monnelll979l: lKentJll989D 



2/ x 1 ( M 



p b {r')(j 2 r (r')- 



dr . 



(15) 

where ybbif) is the bulge surface brightness profile given in 
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Figure 3. The upper and lower panels show the measured bulge 
velocity dispersion along the bulge major and minor axis, respec- 
tively. The solid curve represents our be st-fit model. The filled 
circles represent data of Kormcndv 1 1988) along PA=55° (upper 
panel) and along PA=160° (lower panel). Also plotted are the 
major and min or axis velocity dispersion d ata of IMcElrovl Il983t) 
as smoothed by l^i^row^Perrett^^^u^iJJ^OQ^liopen squares), 
the major axis data of iKoTmend^ &^BendgTJlQggt) (open trian- 
gles), and the major axisdato^f TKonu^n^^CflS^T along PA=38° 



(open circles). The upturn in the curves at 
the influence of the central BH. 



0.01 kpc is due to 



tlBlK Since the projected velocity dispersion is luminosity- 
weighted and since our model isophotes are circular whereas 
the actual bulge isophotes are elliptical, it would seem that 
we perhaps ought to rescale our profiles so that along the 
major axis, cr p (r) is stretched out by ~10% and along the 
minor axis, cr p (r) compressed by the same amount. However, 
we omit this step, since u p {r) is so flat that this rescaling 
has a negligible effect. 

After comparing our results to those obtained with more 
sophisticated assumptions about the bulge geometry and ve- 
locity distribution, we find that the three assumptions listed 
above, especially assumptions (b) and (c), do not appear 
to affect our res ults significantly. We base t his as serti on on 
the discussion in ISimien. Pellet, fc Moniietl (119791) and lKentl 
(1989). The former, for example, demonstrate that impact 
of rotation on the determination of <J p {r) is negligible for 
r < 0.3 kpc. Between 0.3 and 1.0 kpc, the bulge rotation 
profile rises from V c ,b ~ 0.3<t p to 0.5a p . It might be thought 
that this increases the support from rotation and induces 
the slight discrepancy between our models and the data 
at ~ 1 kpc. As reasonable as this sounds, the argument is 
not supported by the results of lKend il989T) . In constructing 
his mass model, Kent assumed, like us, an isotropic bulge 
velocity dispersion but, unlike us, he also used an oblate 
spheroidal model for the mass distribution of the M31 bulge, 
explicitly took into account the axisymmetric nature of the 
disk potential, and allowed for bulge rotation. Kent's veloc- 
ity dispersion profiles are, however, similar to ours both in 



Figure 4. A comparison of the mass profiles for different mod- 
els: our formal best-fit model (solid red line), our best-fit models 
for (M/Lfl) [ j=3.3 (dotted blue line) and the maximum disk (dot- 
long dashed magenta line), the two a xisymmetric models of Ibata 
et al. based on the mass models of iKlvoin. Zhao, fc Somervilld 
(2002) (the dashed and dot-dashed green l ines) , as well as the 
bulge-disk-halo model of lBekki et alJ l200ll) used by Font et al. 
(long da shed cyan line) . The soli d circles represent ma s s estim ates 
made bv lEvans et alJ l200(tl and|E vans fc Wilkinsonl i200(t that 
we used as constraints on our mass models. The first two solid 
circles are mass estimates from planetary nebula and glo bular 
cluster tracers, respectively, from lEvans fc Wilkinsonl l200Cft . The 
last s olid circle is a mass estimate based on M31 satellite galaxies 
from lEvans et all i2000t) . The open squares represent the mass 
estimate based on the observed rotation profile computed using 
the spherical approximation. The open triang le represents an es - 
timate of the mass within 125 kpc derived bv llbata et al.l J2004TI 
while the open circles represent mass estimates from many source s 
which are summarized in Table 6 of [Evans fc Wilkinsor] i2000l) . 
None of the mass estimates shown as open symbols were used as 
constraints on our mass model. 



shape and magnitude. In particular, the Kent model pro- 
files also do not fall off as rapidly as the observed profile. 
To the extent that we have an excelle nt overall ag reement 
with the more rigorous calculations of Kent] lll989l) . we feel 
justified in having made our simplifying assumptions. A clue 
to the origin of the steeper than predicted decline in a p (r) 
may lie in the different ra tes at which the major and minor 
axis profiles fall off. Both ISimien. Pellet, fc MonnetJ dl979tl 
and [Kent] Il988t) note these differences may be indicative of 
velocity anisotropy. 

We note that our neglect of bulge rotation is primarily 
motivated by consistency with our assumption of a spheri- 
cal bulge, and by its small l evel relative to t he bulge veloc- 
ity dispersion. Additionally, IMcElrovl lll983l) noted several 
asymmetries in the rotation curves on the two sides of the 
putative axis of rotation as well as along various different 
position angles. This suggests the determination of a mean 
bulge rotation velocity is not straightforward, which also en- 
courages us to neglect it. 
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3.4 Total Mass Estimates from the Intermediate 
and Outer Halo Regions 

Finally, we consider the constraints on the total mass of the 
system from dynamical studies of globular clusters, plane- 
tary nebulae and M31 satellites. A number of researchers 
have used these tracers to estimate the total mass profile of 
M31 in the intermediate and the outer regions of the halo. 
However, few of these estimates contain a well-motivated 
error bar, which we need for our fitting procedure; many 
authors simply propagate the observational errors through 
their mass estimator and ignore the larger statistical errors 
al together. For the pre s ent pu rp os es, we thus use the r esult s 
of lEvans fc Wilkinson! (I200CT) and lEvans et alJ feOOOl) . who 
estimate the errors using Bayesian statis tics. Using the plan- 
etary nebulae and globular cluster data. lEvans fc Wilkinson! 
(2000) estimate that M(< 31 kpc) = (2.8tj *) x 10 11 M© 
and M(< 40 kpc) = (4.7l^1) x 10 11 M , respectively. Of 
these two, lEvans fc Wilkinson! (|200(J) argue that the mass 
estimate based on globular clusters is more robust because 
these tracers are more uniformly distributed while the plan- 
etary nebulae tend to be concentrated in two regions near 
the o ptical disk. To esti mate the mass further out in the 
halo, lEvans et alJ J2000) analyzed the data for the satel- 
lite galaxies of M31, including the spiral galaxy M33. We 
use their results to estimate the mass at 158 kpc, which 
is the radius of the median object (And II) in the sam- 
ple of 13 M31 satellites, excluding Pegasus and IC 1613. 
These two satellites appear to be outside M31's virial ra- 
dius and thus are not appropriate to include in an analysis 
that assumes statistical equilibrium. We obtain the mass 
estimate by averaging the results of their two velocity dis- 
tributions, and making a small correction from their stated 
total mass to M(< 158 kpc) using their mass profile. This 
gives M(< 158 kpc) = 5.5±^ x 10 11 M©. For each of the 
three halo mass estimates, we convert the given confidence 
intervals to symmetrical error bars in log(M) for use in our 

x 2 fit. 

In Figure 01 we show the above three mass estimates as 
filled circles. These are the values that we use to constrain 
our mass models. For illustrative purposes, however, we also 
plot a number of other mass estimates. The open squares ex- 
tending from ~2 to ~30 kpc trace the mass estimates based 
on the rotational curve data, computed using the spherical 
approximation V? = GM/r. The open triangle at 125 kpc 
designates the mass estimate derived bv llbata et alJ (120041) 
from their preliminary analyses of the dynamics of the gi- 
ant southern stream. We shall consider this further in Pa- 
per II. Finally, the open circles represent mass estimates 
by a number of different r esearchers, taken from Table 6 of 
lEvans fc Wilkinson! i2000T) . 

Figure 0] also shows the results of our best-fit model 
(solid red curve) and two variants, one with the disk mass- 
to-light ratio fixed to (M / Ln)d=3.3 (dotted blue line) and 
the maximum disk model (dot- long dashed magenta line). 
All three curves are in excellent agreement with not only the 
three data points used as constraints, but also in good agree- 
ment with the entire slew of mass estimates shown in the 
figure. The curves, however, do lie systematically below the 
mass estimates computed from the observed rotation profile 
using the spherical approximation. Since the rotation curve 
associated with these models are in excellent agreement with 
the observed profile, this comparison shows clearly the effect 



of neglecting disk flattening. The virial mass (JV/200) corre- 
sponding to the three axisymmetric models spans the range 
6.8— 9.4x10" M© (see Table 0. 



4 PROPERTIES OF ALLOWED MODELS 

The parameters corresponding to our formal best-fit mass 
model for M31 are listed in Table [5] The table also shows 
parameters corresponding to two alternatives to our formal 
best-fit model, one where the disk mass-to-light is fixed to 
(M I Lu)d=S.3 and another where the disk mass is fixed to 
the maximum allowed value. We discuss below the reasons 
for considering these alternative solutions. 

4.1 Allowed Region of Parameter Space 

In general, attempts to model the mass distribution in galax- 
ies are hampered by degeneracies between the various pa- 
rameters. In our specific case, we had expected to encounter 
a linked bulge-disk-halo degeneracy where the halo and the 
disk played off against each other to match the M31 rota- 
tion profile, while the disk and the bulge did the same with 
respect to the inner region dynamics. In fact, our combined 
use of photometric and kinematic data sets greatly reduces 
the available degrees of freedom. 

Of the 8 parameters, two — the disk and bulge mass- 
to-light ratios — are not really independent variables. They 
appear in Equations 1121 and I13H as multiplicative factors 
relating disk and bulge mass profiles to the M31 light pro- 
files. Consequently, once the disk parameters (£0, Rd) are 
determined, the M31 surface brightness profiles fix the disk 
mass-to-light ratio to high accuracy. The same is true of the 
bulge. The light profiles also determine the bulge and disk 
scale radii fairly precisely to 0.61±0.10 and 5.40±0.13 kpc, 
respectively. We therefore fix rj, and Rd to these values after 
doing the surface brightness fit, removing these particular 
degrees of freedom entirely. To esti mate the errors, we have 
adopted the approach advocated bv lPress et al.l (1200111 and 
rescaled the error estimate associated with individual data 
points such that the effective reduced x 2 is unity. As noted 
in H3.ll our fit to the M31's light profile is poor in a formal 
sense using the original photometric errors (the correspond- 
ing reduced \ 2 ls 7) because we do not attempt to fit the 
various visual features present therein. By increasing the er- 
ror bars, we are effectively rendering the physical variations 
in the light curves due to dust lanes, disk warps and spiral 
arms statistically insignificant. We then use this rescaled x 2 
to approximate the uncertainties in the bulge and disk scale 
radii estimates. 

Additionally, we find that the two halo parameters 5 C 
and rh are coupled so strongly that if all other parameters 
are already fixed, setting one effectively sets the other. The 
main reason for this is the small error bars on the rotation 
curve points, which put tight constraints on the halo con- 
tribution there. The constraints on the outer halo force are 
not as strong, as can be seen in Figure 0] 

Moreover, the solutions with the lowest x 2 lie on a thin 
plane corresponding to a specific set of bulge masses Mb. 
The values of this parameter are remarkably robust in that 
forcing large changes in the values of the other parameters 
result in changes of < 10%. This was previously noted by 
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Table 2. M31 Mass Model Parameters for Best-fit, disk with (M/L R )d=3.3, and Maximum disk cases 



Parameter 


Symbol 




Best-fit Model 


<M/L R ) d =3.3 

\ / xx / ft 


Maximum Disk 


Black Hole Mass 


M. 


(10 7 M Q ) 


5.6 


5.6 


5.6 


Total Bulge Mass 


M b 


(10 10 Mq) 


3.3 


3.2 


3.3 


Total Disk Mass 


Ai rf 


(10 10 Mm) 


8.4 


7.2 


13.7 


Total Mass inside 125 kpc 


M(< 125 kpc) 


(10 11 Mq) 


5.6 


6.2 


7.3 


Virial Mass 


M Zoo 


(10 11 Mq) 


6.8 


7.1 


9.4 




M 100 


/1 nil l\/r \ 
(10 Mq) 


/ .O 


1 . ; 


1 n c 
lU.o 


Bulge Scale Radius 


Rb 


(kpc) 


0.61 


0.61 


0.61 


Disk Scale Radius 


w , 
tid 


\ K P C ) 


■ ). 1 


■j . 1 


■J . 1 


Halo Scale Radius 


Th 


(kpc) 


8.18 


7.63 


28.73 


Virial Radius 


^200 


(kpc) 


180.0 


182.3 


200.3 




-^100 


(kpc) 


231.8 


234.8 


262.4 


Halo Density Parameter 


<5c 


(10 4 ) 


27.0 


34.4 


1.54 


Halo Concentration Parameter 


C200 = R200 /rh 




22.0 


23.9 


7.0 


Bulge M/L 


(M/L R )l* 




3.9 


3.9 


4.0 


Disk M/L 


(M/L R )* d * 




3.9 


3.3 


6.3 


Disk Central Surface Density 


So 


(10 s M /kpc 2 ) 


4.6 


3.9 


7.5 


Maximum Rotation Velocity 




km/s 


256.5 


255.9 


256.5 


Fraction of "galactic" baryons 


n m (M b + M d )/(n b M 200 ) 


0.99 


0.86 


1.0 



* We define M/± as the mass enclosed with the sphere of radius such that the mean density inside is Ap c , where A = 100 or 200 

and p c = 277.72/i 2 Mq/ kpc 2 is the present-day critical density. 
** The quoted M/L ratios are based on luminosities that have not been corrected for internal or foreground extinction. 



IWidrow. Perrett. fc Suvul i2003l) . who described the neigh- 
bourhood of their best-fit solution as a "trough [in the M<j- 
M b plane] running parallel to the Md axis." 

This then leaves just two parameters, £0 and one of the 
two halo parameters, to span most of the available parameter 
space. In Figure |5] we show the location of our formal best- 
fit solution and the corresponding \ 2 contours in the So-r?, 
plane. Here we show contour levels of A\ 2 = X 2- Xmin — 1-0, 
2.7, and 4.0, corresponding to 68%, 90% and 95% confidence 
intervals in the case of one parame ter, and 39%, 74% , and 
86% in the case of two parameters jPress et alj EoOl). (Ini- 
tially we regard both parameters as free, but later we will 
assume a value for the disk M/L, resulting in just one free 
parameter.) This plot shows that there is a huge area of the 
plane that produces acceptable fits, although the transverse 
width of this allowed region in the remaining dimensions 
of parameter space is rather small. The importance of the 
formal best-fit solution is thus diluted by this large allowed 
volume, and it is worth seeing whether we can put any more 
constraints on the family of solutions as a whole. 

In the upper panel, we superpose on the contour map 
lines of constant total mass within 125 kpc. The shape of 
the contours shows that the main degeneracy is the tradeoff 
between the disk and halo components; this is not very sur- 
prising, as this tradeoff is a long-standing p roblem in fitting 
galaxy rotation curves iDutton et alj|2005tl . The allowed re- 
gion extends from "maximum disk" models at the upper 
right, to models where the disk is a minor contributor to 
the rotation curve at the lower left. The extension in the 
other direction, across the mass contours, is due to the large 
uncertainty in the outer halo mass, which allows for a large 
tradeoff between and S c . 

The best-fit solution lies on the line for M(< 125 kpc) = 



6.1 x 10 11 M . As already noted. Ilbata et alJ J2004I) have 
carried out a preliminary analysis of the dynamics of the 
giant southern stream, a separate constraint to those exam- 
ined here, and their results require that M(< 125 kpc) > 
5.4 x 10 11 Mq. This constraint reduces the region of allowed 
parameter space by excluding the lower right corner region 
of the A% 2 = 4.0 confidence region, along with the entire 
lower right quarter of the Eo-r^ plane. 

In the middle panel of Figure |S] we superpose lines 
of constant halo concentration C200 on the contour map. 
(Our notation emphasizes we use the radius enclosing a 
mean density of 200 times the critical density to de- 
fine the virial radius and concentration, in contrast to 
some other authors.) Numerical simulation studies suggest 
that in the absence of baryons, the concentration param- 
eter s of ~ 10 12 Mq halos are in the range C200 ~ 7- 
14 (Bullock et alJ l200lt IWechsler et~ail 120021: iDolae et al] 
2004). If anything, the inclusion of baryons ought to make 
the halos more concentrated due to gravitational compres- 
sion of the inner regions of the halo, following the ex- 
pected cooling and pooling of the baryons at the centre 
of the halos fel umenthal, Faber. Flores. fc Primackl Il986l : 
iKlvpin. Zhao, fc Somervillel2002|) 7 On the other hand, heat- 
ing of t he dark matter b y gas clumps during galaxy for- 
mation ( |M^^^Mad|2004l) and by a bar in the central disk 
( Wei nberg fc Kata 120021) may partially counteract this ef- 
fect. For these reasons, we use C200 > 7 as a lower limit 
on the effective dark halo concentration. This condition ex- 
cludes the upper third of the Erj— r% plane, including some 
of the A\ 2 = 4.0 confidence region. 

Finally, the lower panel superposes the contours of the 
fraction of baryons within the virial radius -R200 that are 
"galactic", i.e., within the bulge or the disk component: 
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Figure 5. The cross in the above three panels marks the location 
of the formal best-fit solution (see Table 0, and the solid black 
curves show the corresponding A\ 2 =1.0, 2.7, 4.0 contours. The 
upper panel superposes the contours corresponding to constant 
total mass within 125 kpc (red dashed lines); the middle panel 
shows lines of constant halo concentration C200 (magenta dashed 
lines); and the lower panel shows lines of constant baryon fraction 
(purple dashed lines) as defined in the text. 




0.2 0.4 0.6 0.8 1.0 

Disk surface density (10 9 M e kpc" 2 ) 



Figure 6. The above plot again shows our formal best-fit model 
(cross) and the Ax 2 =1.0, 2.7, 4.0 contours (black solid lines). 
However, we have also combined, for clarity, the three critical 
demarcations from the panels in Figure 151 (see text). We assert 
that the allowed region for physically plausible solutions is the 
section of the Ax 2 <4.0 region that is bounded at the top by 
C200 = 7 (magenta dashed curve) and to the right by / ga i = 1 
(purple dotted line). Additional lines on the plot correspond to 
M{< 125 kpc) = 5.4 X 10 11 M (red dot-dashed curve) and disk 
mass-to-light ratio (uncorrected for extinction) of (M/Ln)^ = 
3.3 (solid blue vertical line). Our constrained, best-fit, physically 
plausible solution is shown as an open circle on the (M/Ln)^ = 
3.3 line. Our constrained, maximum-disk solution is shown by the 
green diamond. 



/g a i = Q m (M d + M b ) / (Q b M 2 oo) ■ In computing this quantity, 
we assume that the bulge and the disk are purely baryonic. 
This is a long-debated assumption even in the case of the 
Milky Way, but the disk and bulge are probably not domi- 
nated by dark matter so that the error from this assumption 
is probably small. Since / ga i must be less than unity, much 
of the lower right half of the En-r/i plane is excluded. For 
Eo > 4.5 x 10 8 M© /kpc 2 , this condition is stricter than that 
of M(< 125 kpc). Figure [S] shows all three constraints on a 
single plot. These constraints, in combination with the A\ 2 
contours, severely restrict the region of allowed solutions. 

The uncorrected disk mass-to-light ratio of our formal 
best-fit solution is (M/Ln)d = 3.9, which is a bit on the high 
side. The mean lumi nosity-weighted B — R color o f M31's 
disk is B - R » 1.60 JWalterbos fc Kennicutlll987l) . Using 
the color-M/ L relations of lBell fc de Jond 11200111) , this corre- 
sponds to a mass-to-light ratio of M/Lr « 3.5. Admittedly, 
the above B — R color is uncorrecte d for reddening due to 
dust while the lBell fc de Jond" j200ll) relations are for dust- 
free quantities. However as they note, the very dust that 
is responsible for the reddening also extinguishes the light 
from the galaxy, increasing the M/L ratio, and to first order 
the two effects cancel out. A more careful analysis suggests 
that the M/L derived using uncorrected colors can be an 
overestimate by as much as ~0.1 dex. We therefore adopt a 
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value of M/Lr ~ 3.3 (uncorrected) for the M31 disk, though 
we recognize that this ratio could be as low as 2.8. 

We can also estimate the corrected or true mass-to- 
light ratio for the M31 disk. The constant foreground red- 
dening inside the Mil ky Way in the direction towards M31 
is E(B - V) = 0.08 jBurstein fe Heiledll98lh whereas the 
internal reddenin g in the disk of M31 is estimated to be 
E(B - V) « 0.25 JWilliams fc HodgdlpOOlTl for a total red - 
dening of E(B - V) « 0.33. Moreover. lBarmbv et aD (l200(t) 
find that for the M31 disk (see their Table 6), E(B - R) w 
1.6x_E(_B — V) « 0.53. Correcting the observ ed B—R color of 
M31 d isk by this amount and then, using the lBell fc de Jonel 
i200ll) correlations yields M/Lr ~1.3. Comparing this to 
the uncorrected M/Lr ss 3.3 implies that the extinct ion in 
R is ~ 1.0 mag, which is in good agreement with the lKen3 
( 1989) estimate of ~0.99 mag for the total extinction. This 
agreement gives us confidence in our adopted value for the 
uncorrected mass-to-light of the M31 disk. 

The relationship between (M/Lr),! and the disk cen- 
tral surface density in our models is (M/Lr)<i = 8.46So,9 
where So, 9 is So in units of 10 9 M0/kpc 2 . In Figure [(J 
we show the vertical line in the So— r~h plane corresponding 
to M/Lr — 3.3. We assert that suitable M31 mass models 
lie on the segment of this line that falls within the existing 
A\ 2 = 4.0 confidence region. In Tabled we give the param- 
eters for the best fit-mass model subject to the additional 
constraint that (M/L_r),j = 3.3. We shall refer to this fit 
as the "M/L-constrained" solution, and we shall use it in 
preference to the formal best-fit solution for the remaining 
discussion. 

Another interesting variant is denoted by the diamond 
located at the intersection of the C200 = 7 and / ga i = 1 con- 
straints. We refer to this as the "maximum disk" solution 
and its properties are also listed in Table [5] This solution 
differs from the traditional maximum disk models in that 
it is required to be consistent with the dynamical data as 
well as the physical constraints on the halo concentration 
parameter and baryonic fraction. Due to the latter two re- 
strictions, the disk in our "maximum disk" solution does not 
dominate the gravitational potential to the same extent as 
in the traditional "maximum disk" models. The disk mass- 
to-light ratio of our "maximum disk" model is M/Lr — 6.3. 

Comparing these two constrained solutions to our un- 
constrained results, we find that the bulge properties, such 
as the model mass as well as the realization of the bulge ve- 
locity dispersion and surface brightness profiles, are nearly 
identical. The altered surface density in the disk is offset 
by the halo becoming denser and more concentrated for the 
constrained (M/Lr),i — 3.3 model, and the reverse for the 
maximum disk model, as seen in Figure [5] and Table [5] The 
overall effect is that the total rotation profiles for the con- 
strained and the unconstrained cases are nearly the same. 
In the M/L-constrained model, the halo contribution to the 
potential is greater th an that of the disk at all radii, a prop- 
erty that according t ol Widrow. Perrett. fc Suvul i2003t) and 
IWidrow fc Dubinskil (I2OO5F is likely to delay or even sup- 
press the onset of bar formation. In contrast, the gravita- 
tional contribution of the disk in our maximum disk model 
is greater than that of the other two components over the 
radial range 5 kpc < r < 25 kpc, and based on their N-body 
simulation studies. lWidrow fc Dubinskil j2005T) indicate that 
such models are highly susceptible to the formation of strong 
persistent bars. 



4.2 Comparing to Other Mass Models 

iKend i 19891) was one of the first to construct a detailed 
mass model for M31 using both photometric and kinematic 
data. Comparing our M/L-constrained model (Tabic |2J to 
his small bulge model, we find that his bulge properties are 
in excellent agreement with ours: Kent's bulge is a factor 
of ~ 1.2 more massive than ours and the corresponding 
mass-to-light ratio is (M/LR)b = (3.7 ± 0.4) , compared 
to our value of 3.9. Kent's bulge is represented as an oblate 
spheroid whereas our bulge is spherically symmetric, while 
other differences in our assumptions are detailed in § 3.1 
above. In view of these differences in the models, the level 
of agreement in our results for the M31 bulge supports our 
prior assertions regarding the robustness of the bulge pa- 
rameters. Kent's disk, on the other hand, is approximately 
twice as massive as that of our M/L-constrained model, and 
the corresponding mass-to-light ratio of (7.6ZIZ0.8) 1 is also 
larger than our value of 3.3. As indicated by his Figure 2, 
the disk is the primary contributor to the rotation profile 
over the range 5 kpc < r < 25 kpc whereas in our preferred 
model, the disk does not dominate at any radii. However, the 
bulge and disk properties, including the mass-to- light ratios, 
of Kent's model are similar to those of our maximum-disk 
mode l . Am ong the differences in our two schemes are: (a) 
IKend ( 1989) fixed the central surface brightness of his disks 
to values lower than what we find; and (b) he adopted a con- 
stant density dark matter halo whereas we assume an NFW 
functional form. The differences in our preferred results seem 
primarily due to assumption (b). Given these different as- 
sumptions underlying our two approaches, and the disk-halo 
degeneracy discussed in § 4.1, the disagreement in our disk 
resu lts, at a factor ~ 2 level, i s not s urprising. 

IWidrow. Perrett. fc Suvul J2003ll developed their M31 
model much like we did; they sought out best-fit solu- 
tions by minimizing a composite % 2 statistic, and some 
of the datasets we use are taken from their paper. How- 
ever, their approach differs from ours in that they mod- 
eled the disk, bulge and halo components of the galaxy 
as distribution functions, which cannot be specified ana- 
lytically. In practise, this meant that their halo was rep- 
resen t ed by a low ered Evans model llKuiiken fc Dubinskil 
I1994 lEvans! Il993h instead of an NFW-like profile, their 
bulge by a lowered isothermal sphe re or a King model 
<Kindll96fi iBinnev fc Tremainelll987ll . and their disk by a 
lKTrilken"fc]3ulTmskin^)94l) model. This resembles our disk 
in that its surface density falls off exponentially with radius; 
unlike us, they also take into account the finite thickness of 
the disk. Widrow et al. list the properties of their best-fit 
model as well as a number of variants. In terms of its charac- 
teristics for r < 30 kpc, our best-fit M/L-constrained model 
is similar to their preferred model (Model A), as can be 
seen by comparing their Figure 4 to our Figure [5] There 
are, however, some quantitative differences. Their decompo- 
sition of the surface brightness into bulge and disk compo- 
nents gives a bulge-to-disk light ratio of 0.58, versus our ratio 
of 0.38. This may stem from their different functional forms 

1 The se values are slightly different from those quoted in iKentj 
( 1989) because they have been rescaled to account for the 
differences in our assumed distances to M31, and converted 
from r-band to R-band using M31 colors and (r — R)q = 0.18 
( http:/ /www. sdss.org/dr4/algorithms/sdssUB VRITransform.html I 
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for the disk and bulge, or it may b e because they use only 
the gl obal surface brightness from IWalterbos fc Kennicutd 
( 1987), which does not constrain the bulge region very well. 
Their bulge-to-disk mass ratio, in contrast, is slightly smaller 
than ours: 0.37 versus 0.44. As a result, their bulge M/L 
is 2.7, 30% smaller than ours, and their disk M/L is 4.4, 
33% larger than ours. The contrast between their bulge and 
disk M/L is not very reasonable, given that the observed 
disk color is slightly bluer than that of the bulge, unless 
there is substantial dark matter in the disk. As for the halo, 
Widrow et al. replace their lowered Evans model after the 
fact with an NFW halo that closely matches the inner struc- 
ture of their original halo. This NFW halo has C200 = 11.5, 
r h = 19.5 kpc, R200 = 224 kpc, and M 20 o = 1.3 x 10 12 M Q . 
The resulting total mass of M31 is roughly a factor of 1.8 
greater than that of our system, and is somewhat on the 
high side in light of the mass constraints displayed in Fig- 
ure 0] Nevertheless, the disk and halo properties still put 
their model within the allowed region in Figure |SJ and in 
general the differences between their model and ours are 
not very substantial given the large differences between our 

fun ctional forms. 

iKlypin. Zhao, fc Somervilld <|2002T) also constructed a 
mass and light model of M31, in order to test the consis- 
tency of large disk galaxies such as the Milky Way and M31 
with the cuspy halos of the LCDM paradigm. An interesting 
featu re of this model is the inclusion of adiabatic contrac- 
tion ( Blu menthal. Faber. Flores. fc Primacklll98rJ) of their 
initially NFW halos by the gravity of the disk and the bulge. 
To fit the baryonic component and the halo, they use a 4- 
component model containing a flattened bulge, a bar, an 
exponential disk, and the halo. This contains about 14 pa- 
rameters in total. They do not actually state the parameter 
values or how the parameters were chosen, but they assert 
that these achieve a good match to the data without the 
need for careful optimization. Their treatment of the bulge 
and disk dynamics is even more approximate than ours; the 
central mass profile is calculated from an isothermal sphere 
approximation, and the disk rotation curve apparently does 
not take flattening into account. For this reason, some differ- 
ences from our results are to be expected. Nevertheless, we 
find a remarkable level of agreement between their Model 
CI and our best-fit constrained model, as can be seen by 
comparing their Figure 4 to our Figure [5] One primary dif- 
ference is that their bulge-disk decomposition is different, 
with the mass ratio between the two only 0.27 whereas in 
our model it is 0.44. This could be a product of the dy- 
namics modeling or of the different bulge profiles. A second 
difference is that, correcting their unextincted luminosity to 
extincted using the transformations they suggest, they ob- 
tain a total galaxy R luminosity of 3.5 x 1O 1O L0 or 16% larger 
than ours; we do not understand the origin of this difference, 
since we are using nearly the same data and they have not 
raised the |Walterbos fc Kennicuttl il987f) surface brightness 
by 0.1 mag as we have. Together with their total baryonic 
component mass, which is smaller than ours by about 10%, 
the higher luminosity implies a total galactic (bulge+disk) 
mass-to-light ratio of M/Lr = 2.5 (extincted) whereas in 
our model it is 3.5. A third difference is that they use a 
somewhat more massive halo. Overall, we estimate that the 
mass of their total system within 125 kpc is approximately 
0.9-1.0 x 10 12 M Q . In our best-fit constrained model, the 
total mass within 125 kpc is roughly a factor of 1.5 less. 



This difference is probably within the large errors on the 
halo mass and is not necessarily a product of our different 
modeling technique. 

iBermanl l|200ll) and iBerman fc Loinardl <2002l) mod- 
eled the gravitational potential in the inner ~ 10 kpc of 
M31, using a rotating triaxial bulge. By comparing to the 
molecular gas flows seen in high-resolution CO velocity 
maps, they found a best-fit bulge mass of 2.3 x 10 10 M© 
within 3.5 kpc (see Ber man fc Loinardl 120021) . In our three 
models in Table |5J the bulge mass w ithin this radius is 
2.4 x 10 10 M Q . The agreement with the lBerman fc Loinardl 
(2002) results is excellent, especially considering the com- 
plete independence of the datasets involved, and the sim- 
plifying assumptions we have made in the treatment of the 
bulge which (as em phasized above) is n ot the principal con- 
cern o f our models. IBermanl i200l]) and lBerman fc Loinardl 
( 2002) also estimate the disk and the total mass within 3.5 
kpc to be 1.2 x 10 10 M Q and 3.7 x 10 10 M , respectively. 
Given that the latter is nearly the same as t he co mbined 
disk+bulge contribution, IBerman fc Loinardl <|2002T) argue 
that the gravitational contribution of the halo within 3.5 
kpc must be negligible. Our analysis, however, does not 
support this latter conclu sion. We find that although the 
IBerman fc Loinardl <l2002t) disk mass estimate is compara- 
ble to that in our formal best-fit and M/L-constrained mod- 
els (1.1 and 1.0 x 10 10 Mq, respectively), we require a total 
mass within 3.5 kpc of 5 x 10 10 Mq, a factor of 1.35 times 
larger, to account for the rotation curve data. Consequently, 
the dark halo contributes roughly a third of the total gravity 
within 3.5 kpc. 

In th eir an alyses of the stream dynamics, 
Ibata et al.| <| 2004|) consid e red two variants of the 
Klypin. Zhao, fc Somervillel ]2002l) Model CI in which 
the bulge and disk properties were kept fixed to the values 
specified by Klypin et al., but the halo properties were 
allowed to vary with the aim of minimizing discrepancies 
between the orbit of the progenitor and that of the giant 
stream. Dependin g on which stream constraints they used, 
llbata et alJ i2004T) iterated to two models. The rotational 
velocity and the total mass profiles (computed from the 
rotation curve using the spherical approximation) for these 
two models are shown as short dashed and dot-dashed 
(green) lines in lower right panel of Figurc|5|and in Figure|l| 
As illustrated in Figure |21 the rotatio n curve for the mo del 
dubbed "high mass" (short das hed) bvllbata et al] (120041) is 
in onflict with the observations, llbata et al.fil2004T) note this 
and reject the model. The characteristics of the alternate 
"low mass" model (dot-dashed curve) are broadly similar to 
those of our constrained best-fit axisymmetric case, and for 
the purposes of orbit calculations, we expect that the two 
will give br o adly s imilar results. However, the fact that the 
llbata et al] i2004T) system has 20% more mass within 125 
kpc will, for a given set of initial conditions, lead to small 
differences in the apocenter. We discuss this point further 
in the next section. 

iFont et "all (120041) . in their recent study, a dopt the ax- 
isymm etric bulge-disk-halo M31 mass model of lBekki et al] 
( 2001) with potential: 



*(r) =$b(r) + Mr)+$h{r), (16) 
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where, 



Mr) 
Sd(r) 



GM b 
n + r' 



GMd 



^R 2 + (a + Vz 2 + b 2 ) 2 ' 
\v h 2 \n{r 2 + r 2 ), 



(17) 
(18) 

(19) 



r is the distance from the center of the galaxy, 7? and z are 
polar coordinates aligned with M31's disk plane, r c — 12 
kpc, V h = 186 km s"\ M d = 1.3 x 10 11 M , a = 6.5 kpc, 
b = 0.26 kpc, M b = 9.2 x 10 10 M , and r b = 0.7 kpc. 
For comparative purposes, the rotation curve and the total 
mass profiles (computed from the rotation curve using the 
spherical approximation) for this model are plotted in the 
lower right panel of Figure [5] and in Figure 2] as long-dashed 
(cyan) curves. The model rotation curve is systematically 
higher than the observed profile, with the two being least 
discrepant in the outer disk region where the model curve 
is approximately 20 km s" 1 higher, and most discrepant in 
the inner disk. Th e latter is primarily due to the large bulge 
mass assumed bv lBekki et alJ (1200 J) . Their bulge mass is 
a factor ~ 2.5 larger than ours. Additionally, the Bekki et 
al. model overestimates the the system mass beyond ~ 80 
kpc. On the whole, the Bekki et al. model describes a galaxy 
with a much deeper potential than that indicated by the 
observations or that described by our preferred M31 model. 
As we shall illustrate in the next section, this has significant 
implications for orbit calculations. 

To summarize, the models of the M31 potential pre- 
sented in this paper seem to be more consistent with the 
observations than several others presented in the literature, 
and generally consistent with others of greater complexity 
while being easier to use for calculating orbits. There is, 
however, still a great deal of uncertainty in the actual poten- 
tial. Our model, for example, is more strongly constrained in 
some ways by consideration of the disk mass-to- light ratio, 
the halo concentration, and the galaxy baryonic fraction, 
than by the actual observations. This uncertainty has bear- 
ing on studies of the response of the disk to perturbations, 
the response of dark halos to the assembly of galaxies, and 
the orbits of satellite galaxies, among other issues. We expect 
that the model can be improved in the future by the addi- 
tion of more halo tracers, s uch as the halo stars detected by 
iGuhathakurta et alJ l)2005h and Chapman et al. (in prepa- 
ration); by using coherent stellar streams, such as the giant 
southern stream and the s tream in the vicinity of NGC 205 
l|McConnachie et alJl2004Tl : by the use of deep infrared sur- 
face brightness measurements, which are less sensitive to 
dust and stellar population variations than the optical data 
used here; by the incorporation of non-axisymmetric fea- 
tures such as the oblong/triaxial bulge and the spiral fea- 
tures and warping of the disk; and by properly treating the 
finite thickness of the disk. 



5 SAMPLE ORBITS IN SPHERICAL AND 
FLATTENED POTENTIALS 

With our mass models defined, we are now in the position 
to study the orbits giving rise to the giant southern stream 
in M31. The results of our study are described in detail in 



Paper II llFardal et alJl2005l) . Here, we merely present a few 
examples that illustrate the sensitivity of the orbits to the 
mass models used. 

To facilitate comparison, we start by adopting initial 
conditions similar to those used bv lFont et alJ i2004fl These 
authors compute an orbit that is in rough agreement with 
the kinematic observa t ions o f the giant southern stream 
bv iMc Connach ie et alJ (I2003T) . lLaw. Johnston. & Maiewski 
l|2005h . llbata et alJ l|2004l 'l. and IGuhathakurta et alJ ^OolT 
We assume here, with these authors, that the stream and 
the progenitor orbit are coincident. (This condition will be 
modified in Paper II). Font et al. define a coordinate sys- 
tem (x, y, z) such that the x-y plane is in the sky, with x 
increasing to the E and y to the N, and z is parallel to the 
line of sight. Font et al. locate the starting point of their 
orbit in the centre of the southern stream. For a detailed 
discussion of the stream geometry a nd of the choice o f ini- 
tial conditions, we refer the reader to lFont et al] i2004f) and 
our Paper II. 

We slightly optimize the Font et al. initial conditions 
to better match the observations of the stream in our M/L- 
constrained disk potential. Specifically, we set 
xo = 14.2 kpc, Vxo — — 50 km s~\ 
yo = —33.2 kpc, v y o = 86.02 km s _1 , 
Zq = 60 kpc, v z q = — 158 kms" 1 . 

For the comparisons here, the details of this optimiza- 
tion procedure are not important. Issues surrounding the 
comparison to observations will be discussed in Paper II. 
The primary thing to note is that the apocenter of the orbit 
must occur approximately at Field 1, since the z velocity of 
the stream relative to M31 is nearly zero there. 

The main difference between the initial conditions that 
we use here and those used by Font et al. is that our ini- 
tial velocity is not as high as theirs, because our potential 
for M31 is not as deep as the one they use. As a simple 
illustration of this point, Figure |7] shows test particle or- 
bits calculated forward and backward from our new initial 
conditions, but computed in two different potentials: our 
M/L-constrained mass model (solid curve) and the simple 
singular isothermal sphere halo (dashed curve) with poten- 
tial given by 



$(r) = V 2 In (r), 



(20) 



where V c = 200 km s -1 is the circular velocity. The latter 
is a good approximation to the model used by Font et al. in 
the outer halo. The key thing to note in Figure |7] is that in 
the deeper singular isothermal potential, the test particle is 
unable to climb out as far before reaching apocenter when 
its orbit is traced backward in time, and it falls short of 
Field 1 by approximately one degree. To guarantee that the 
apocenter of the test particle's orbit falls in Field 1 requires 
endowing the test particle with higher initial velocities. Since 
the radial velocities along the stream are observed quanti- 
ties, the differences in the initial velocities will eventually 
result in a different assumed tilt of the orbital plane, and 
thus different projections of subsequent lobes of the orbit. 

In Figure |H| we examine the effect of the strength of 
the flattened disk potential. The dotted line shows the same 
results for the M/L-constrained model as in Figure □ The 
solid line shows the results for the spherical disk mass model 
with M/L-constrained model parameters. (We remind the 
reader that the halo and the bulge components of these two 
models are exactly the same, and even the radial mass pro- 



© 2006 RAS, MNRAS 000.HHT71 



Investigating the Andromeda Stream: I 15 



M/L-constrained 
Isothermal halo 





M/Z.-constrained 
Isothermal halo 




6 4 2 -2 -4 -6 

EW position i; (degrees) 
Figure 7. Comparison of test particle orbits with exactly the 
same initial conditions but in two different potentials: the simple 
singular isothermal sphere potential discussed in the text (dashed 
curve) and the potential corresponding to our M / //-constrained 
model (solid curve). The orbits are started in the centre of the 
southern stream with initial velocity such that the orbit in our 
potential turns around near Field 1. Red and blue curves rep- 
resent orbit integrations forward and backward in time, respec- 
tively. (Top panel): shows the orbit in sky coordinates, § and r/. 
(Middle panel): plots the line-of-sight depth (in kpc) of the orbit 
versus the £ position in the sky. (Bottom panel): plots the varia- 
tion of radial velocity (in km s ) with £ position in the sky. The 
symbols in the three panels are as follows: star (M31), triangle 
(M32), plus (NGC 205), cross (A nd VIII), diamond (Field a3), 
squares iMcConnachie et al 1l2003l stream fields). 
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files of the disk component ALj(< r) for the two models are 
the same. The only difference is that in the spherical disk 
case, the disk mass is distributed in a spherically symmet- 
ric fashion, as opposed to planar and axisymmetric. Con- 
sequently, the M/L-constrained spherical disk model is, as 
the name suggests, a spherically symmetric analog of our 
axisymmetric M/L-constrained model.) The result for the 
spherical disk model is essentially indistinguishable from the 
flattened disk model over the duration of the radial lobe that 
includes the stream. However, the directions of the preceding 
and subsequent lobes are different. These parts of the trajec- 
tories are very sensitive to the differences in the geometry 
of the disk potential at pericenter, since the test particle 
plunges all the way in to where the disk contribution to 
the total potential is significant. The gravitational pull from 
the disk, pulling the test particle towards the plane of the 
disk, results in an orbit that bends more strongly than in 
the spherical case. The dashed line uses the results from 
our "maximum disk" model. Again, the primary effect is a 
change in the direction of the previous and subsequent radial 
loops relative to the orbit in our fiducial M/L-constrained 
potential although in this case, the slightly deeper potential 
also results in a slightly smaller radial excursion at apocen- 
ter. Both of these effects are due to the stronger gravity of 
the disk in the maximum-disk model. 

In summary, the radial mass profile has a significant 
effect on orbits in the halo, and is the primary factor in de- 
termining the length of the stream given the initial velocity, 
or alternatively the velocity gradient given the length of the 
stream. The geometric shape of the mass profile is a less sig- 
nificant factor, at least with the tests performed here, but 
it does significantly affect the direction of subsequent lobes. 
If part of the objective for computing such orbits is to try 
to ascertain where the progenitor might be, then ensuring 
that the orbits are being computed in a realistically shaped 
potential is important. We have neglected the possible flat- 
tening of the halo, both here and in Paper II, but this is 
clearly an interesting direction to explore at some point. 



6 SUMMARY 

The following are several significant results from this work 
that are important to highlight: 

• We have presented a simple analytic model for the mass 
density, light density, and gravitational potential of M31, 
which is easy to use for the purposes of orbit calculations. We 
apply this potential in Paper II to the problem of estimating 
the orbit and other properties of the giant southern stream 
in M31's halo. 

• Our new potential does a better job matching the dy- 
namics of M31 than other simple potentials used recently 
for the dynamics of the stream. It does a comparable job, 
and in fact is fairly similar in its properties, to several more 
sophisticated but non-analytic models in the literature. 

• The eight parameters in our mass and light model span 
a large but low-dimensional region in parameter space. Phys- 
ical constraints on the fraction of halo baryons in the galaxy, 
the halo concentration, and the stellar mass-to-light ratio 
restrict the parameters to a much smaller region. Our pre- 
ferred model uses the M/L ratio of the disk as a further 
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EW position i; (degrees) 

Figure 8. The above three panels show a comparison between 
the models with different strengths of the disk potential. All pan- 
els, symbols, and colours are the same as in Figure [7] however, 
the dashed line now represents the "maximum disk" model, the 
dotted line represents the M/ L-constrained model, and the solid 
line represents the M/L-constrained "spherical disk" model (see 
text for details). All orbits are calculated using the same initial 
conditions as in Figure Q Note that the time intervals for the 
different orbits are not the same. 



constraint. The parameters of this model and its important 
physical quantities are summarized in Table [5] 

• The radial mass profile has a significant effect on orbits 
within the giant southern stream, and is the determining 
factor in the length of the stream for a given initial velocity. 
The geometric shape of the mass profile does not play as 
significant a role in determining the orbit as the radial mass 
profile, at least within the scope of possible variations of 
the disk surface density. However, it does have an effect on 
the direction of subsequent lobes of the orbit. This becomes 
important if one wishes to construct orbits for the purposes 
of locating the progenitor of the stream. 
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